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Abstract 


The  joined- wing  configuration  that  was  published  by  Wolkovich  in  1986  has 
been  studied  by  many  researchers.  The  joined-wing  airplane  is  defined  as  an 
airplane  that  incorporates  tandem  wings  arranged  to  fonn  diamond  shapes  in  both 
top  and  front  views.  The  joined-wing  can  lead  to  increased  aerodynamic 
perfonnances  as  well  as  reduction  of  the  structural  weight.  However,  the  joined- 
wing  has  high  geometric  nonlinearity  under  the  gust  load.  The  gust  load  acts  as  a 
dynamic  load.  Therefore,  nonlinear  dynamic  (transient)  behavior  of  the  joined- 
wing  should  be  considered  in  structural  optimization.  In  previous  researches, 
linear  dynamic  response  optimization  and  nonlinear  static  responses  optimization 
are  performed.  It  is  well  known  that  conventional  nonlinear  dynamic  response 
optimization  is  extremely  expensive.  Therefore,  in  this  research,  nonlinear 
dynamic  response  optimization  of  a  joined-wing  is  carried  out  by  using  equivalent 
static  loads.  The  concept  of  equivalent  static  loads  is  expanded  and  newly 
proposed  for  nonlinear  dynamic  response  optimization.  Equivalent  static  loads  are 
the  load  sets  which  generate  the  same  response  field  in  linear  static  analysis  as  that 
from  nonlinear  dynamic  analysis.  Therefore,  nonlinear  dynamic  response 
optimization  can  be  conducted  by  repeated  use  of  linear  response  optimization. 
For  the  verification  of  efficiency  of  the  proposed  method,  a  simple  nonlinear 
dynamic  response  optimization  problem  is  introduced.  The  problem  is  solved  by 
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using  both  the  equivalent  static  loads  method  and  the  conventional  method  with 
sensitivity  analysis  using  the  finite  difference  method.  The  procedure  for 
nonlinear  dynamic  response  optimization  of  a  joined- wing  using  equivalent  static 
loads  is  explained  and  the  optimum  results  are  discussed. 
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1  Introduction 


The  joined-wing  is  an  innovative  aircraft  configuration.  The  joined-wing  may  be 
defined  as  an  airplane  that  incorporates  tandem  wings  arranged  to  form  diamond  shapes  in 
both  top  and  front  views.  Wolkovich  proposed  a  joined-wing  design  with  potential  weight 
reduction  and  aerodynamic  benefits  as  early  as  1986.  (1)  The  joined-wing  has  the 
advantage  of  a  longer  range  and  loiter  than  those  of  a  conventional  wing.  Generally,  the 
weight  of  the  joined-wing  aircraft  is  lighter  than  that  of  a  conventional  wing.  Fig.  1 
shows  a  general  joined-wing  aircraft  where  the  fore-wing  and  aft-wing  are  joined.  Miura, 
Shyu  and  Wolkovich  employed  an  optimization  method  to  study  the  effects  of  joined-wing 
geometry  parameters  on  structural  weight.  (2)  Gallman,  Smith  and  Kroo  offered  many 
recommendations  for  the  design  methodology  of  a  joined-wing.  13 '  They  used  the  fully 
stressed  design  (FSD)  for  optimization.  Blair  and  Canfield  initiated  nonlinear  exploration 
on  a  joined-wing  configuration  in  2005.  (4)  Air  Force  Research  Laboratories  (AFRL)  have 
been  developing  an  airplane  with  a  joined-wing  to  complete  a  long-endurance  surveillance 
mission.  <4‘8)  Lee  et  al.  performed  linear  dynamic  response  structural  optimization  of  a 
joined-wing  using  equivalent  static  loads.  <9)  They  considered  the  dynamic  effect  of  the 
joined-wing  in  optimization.  In  2007,  Kim  et  al.  performed  nonlinear  static  response 
optimization  of  a  joined-wing  using  equivalent  loads. <10) 

In  previous  researches,  (4'  101  it  is  certain  that  the  joined-wing  has  high  geometric 
nonlinearity  under  the  gust  loading  conditions  due  to  the  specific  shape  of  the  joined-wing. 
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And  structural  optimization  of  the  joined-wing  has  been  performed  from  the  viewpoint  of 
nonlinear  static  response.  However,  real  forces  act  dynamically.  Especially,  the  gust 
loads  are  the  most  important  loading  conditions  when  an  airplane  wing  is  designed.  The 
gust  is  the  movement  of  the  air  in  turbulence  and  the  gust  load  has  a  large  impact  on  the 
airplane.  (1112)  The  gust  loads  generate  various  dynamic  effects  on  the  aircraft  wing. 
Therefore,  the  nonlinear  dynamic  effect  of  the  joined-wing  should  be  considered  in  the 
optimization  process.  However,  it  is  very  difficult  to  optimize  a  joined-wing  considering 
the  nonlinear  dynamic  effect.  The  reason  is  because  the  conventional  optimization 
method  is  not  efficient  for  nonlinear  dynamic  response  structural  optimization. 

The  calculation  of  the  sensitivity  from  nonlinear  dynamic  analysis  is  fairly  difficult. 
This  is  due  to  the  great  number  of  nonlinear  dynamic  analyses  required  for  the  calculation 
of  the  sensitivity.  Therefore,  the  conventional  gradient  based  optimization  method  is  not 
useful  for  nonlinear  dynamic  response  optimization.  (13'15)  The  non-gradient  based 
optimization  method  such  as  the  response  surface  method  can  be  used  for  nonlinear 
dynamic  response  optimization.  However,  the  method  has  several  disadvantages  such  as 
the  limit  of  the  number  of  design  variables  and  inaccuracy  of  the  solution.  <16)  In  this 
research,  the  equivalent  static  loads  (ESL)  method  is  introduced  for  nonlinear  dynamic 
response  optimization.  Until  now,  the  equivalent  static  loads  method  has  been  used  for 
linear  dynamic  response  optimization  and  nonlinear  static  response  optimization.  1 16-21  * 
The  concept  of  ESL  is  expanded  to  nonlinear  dynamic  response  optimization. 

Equivalent  static  loads  are  defined  as  the  linear  static  load  sets  which  generate  the  same 
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response  field  in  linear  static  analysis  as  that  from  nonlinear  dynamic  analysis.  Therefore, 
if  equivalent  static  loads  are  used  as  applied  loads,  the  same  responses  from  nonlinear 
dynamic  analysis  can  be  considered  throughout  linear  static  response  optimization.  It  is 
well  known  that  nonlinear  dynamic  analysis  is  quite  expensive.  On  the  other  hand,  linear 
static  analysis  is  not  costly  and  the  linear  static  response  optimization  theory  is  well 
established.  Equivalent  static  loads  are  made  to  reduce  the  number  of  nonlinear  dynamic 
analyses.  Moreover,  because  the  method  is  gradient  based  optimization,  the  solution  is 
exact.  A  detailed  explanation  will  be  introduced  in  Section  2. 

As  a  small  example,  a  cantilever  plate  is  optimized  under  equivalent  static  loads 
transformed  from  a  dynamic  load.  The  results  are  compared  with  those  of  a  conventional 
method  where  the  finite  difference  method  is  employed  for  sensitivity  calculation.  A 
joined-wing  is  optimized  under  dynamic  gust  loads.  The  gust  loads  are  considered  as 
external  loads  in  nonlinear  dynamic  response  optimization.  The  gust  loads  for  a  joined- 
wing  have  been  calculated  by  the  researchers  of  the  AFRL.  <4)  Static  loads  for  the  gust  can 
be  generated  from  an  aeroelastic  model  which  uses  the  Panel  method.  (IM2)  It  is  difficult 
to  identify  the  exact  dynamic  gust  load  profile.  Therefore,  the  static  gust  loads  from 
AFRL  are  transformed  to  dynamic  loads  using  the  1 -cosine  function.  <U)  ABAQUS  6.7  <24) 
is  employed  for  nonlinear  dynamic  analysis  and  GENESIS  9.0  <25)  is  used  for  linear  static 
optimization. 
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2  Nonlinear  dynamic  structural  optimization  using 


equivalent  static  loads  (NDROESL) 


As  mentioned  earlier,  nonlinear  dynamic  response  structural  optimization  is  quite 
difficult  even  with  the  modem  computer  system.  Nonlinear  analysis  considering  time  is  a 
lot  more  expensive  than  nonlinear  static  analysis.  This  disadvantage  is  fatal  for  structural 
optimization  using  the  gradient  based  optimization  method  because  the  calculation  of 
sensitivity  needs  a  large  number  of  nonlinear  dynamic  analyses.  On  the  other  hand,  the 
approximation  methods  such  as  the  response  surface  method  (RSM)  are  easy  to  use; 
however,  they  have  a  limit  on  the  number  of  design  variables  and  the  solutions  are  not 
exact.  (16)  The  equivalent  static  loads  method  is  a  new  and  efficient  method  that 
overcomes  those  weaknesses.  In  this  section,  the  concept  and  the  calculation  of  the 
NDROESL  method  are  explained. 

2.1  Problem  formulation  of  nonlinear  dynamic  response  structural 
optimization 

The  formulation  for  the  nonlinear  dynamic  response  optimization  can  be  expressed  as 
follows: 


Find 

b  e  Rm 

(la) 

to  minimize 

/(b) 

(lb) 
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subject  to 

M(b)Zjv(0  +  K  v(b,ziV(f))z  v(?)  -  f(f)  =  0 

(lc) 

t  =  1 ,...,« 

gj(b,zN(t))<0,  j  =  l,—,l 

(Id) 

biL^bi  <btu,  i  =  l,—,m 

(le) 

M  is  the  mass  matrix  which  is  the  function  of  the  design  variable  vector  b  .  K  is  the 
stiffness  matrix  which  is  the  function  of  the  design  variable  vector  b  and  the  nodal 
displacement  vector  z ,  and  z  is  the  acceleration  vector.  The  subscript  N  means  that 
the  response  is  obtained  from  nonlinear  analysis.  Eq.  (lc)  is  the  governing  equation  of 
nonlinear  dynamic  analysis  using  the  finite  element  method.  <22' 23)  The  constant  n  is  the 
total  number  of  the  time  steps.  The  constant  /  and  m  are  the  total  number  of  the 
constraints  and  design  variables,  respectively.  f(7)  is  the  external  load  vector  at  the  t  th 
time  step.  biL  and  bjL,  are  the  lower  bound  and  upper  bound  of  the  i  th  design  variable, 
respectively. 

As  mentioned  earlier,  dynamic  response  optimization  has  many  time  dependent 
constraints.  As  shown  in  Eq.  (Id),  the  total  number  of  time  dependent  constraints  is  n  x  / . 
Moreover,  the  calculation  of  the  sensitivity  considering  the  incremental  step  is  extremely 
difficult.  Therefore,  it  is  rare  to  perform  nonlinear  dynamic  response  structural 
optimization  for  large  scale  problems. 

2.2  Calculation  of  the  equivalent  static  loads 

The  equivalent  static  loads  (ESL)  are  defined  as  the  static  loads  which  generate  the  same 
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response  fields  as  those  under  a  dynamic  load  at  an  arbitrary  time  of  dynamic  analysis. 
According  to  the  finite  element  method,  (22‘21)  the  equilibrium  equation  of  a  structure  in  the 
time  domain  with  nonlinearity  is 

M(b)zv(f)  +  K(b,Z7V(0)zw(0  =  f(t) 

(t  =  0,1, 2,...,n) 

zN  at  all  the  time  steps  is  obtained  from  Eq.  (2).  The  equivalent  static  load  for 
displacements  is  defined  as: 

X:q{s)  =  KL{b)iN(t) 

0  =  0,1,2,...,/?) 

where  new  notation  s  is  exactly  matched  with  t  in  Eq.  (2).  The  reason  to  use  the 
notation  s  is  that  Eq.  (3)  is  not  defined  in  a  dynamic  region  but  in  a  static  region.  In 
other  words,  t  =  i  is  equal  to  s  =  i  and  the  total  number  of  5  is  n  .  Therefore,  n 

equivalent  static  loads  are  obtained  from  Eq.  (3).  fj?0)  is  the  equivalent  load  vector  for 

displacement  at  each  time  step,  K L  is  the  linear  stiffness  matrix  and  zN{s)  is  the  nodal 
displacement  vector  from  Eq.  (2).  f lq(s)  is  used  in  Eq.  (4)  which  is  the  equation  of  linear 
static  analysis  as  follows: 

Ki(b)z£(j)  =  f*(j)  (4) 

where  the  nodal  displacement  vector  zL(s )  has  the  same  values  as  the  nonlinear  nodal 
displacement  vector  zN{t)  in  Eq.  (2)  at  an  arbitrary  time.  Therefore,  if  the  equivalent 

static  load  f *g(s)  is  used  as  an  external  load  in  linear  static  response  optimization,  the 

same  displacements  as  the  nonlinear  dynamic  response  can  be  considered  in  linear 
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response  optimization.  The  equivalent  static  loads  are  used  as  multiple  loading  conditions 
for  linear  static  response  optimization  and  Fig.  2  presents  this  process. 

Although  the  load  f *q{s)  can  generate  the  same  displacements  as  the  nonlinear 

displacements  at  all  the  time  steps,  it  does  not  generate  the  same  stress  responses  because 
the  relationships  between  the  strain  and  displacement  as  well  as  the  strain  and  stress  have 
nonlinearity.  Thus,  the  equivalent  static  loads  for  the  stresses  are  separately  calculated. 
The  stress  response  <sN(t)  is  obtained  from  Eq.  (2)  of  nonlinear  dynamic  analysis.  The 
obtained  stresses  are  used  as  initial  stresses  of  linear  static  analysis.  The  equivalent  static 
loads  for  stresses  are  calculated  as  follows: 

Ki(b)z"W  =  -f/(ow(0)  (5) 

f“Cs)  =  K,.(b)z?(s)  (6) 

where  f“  (s)  is  the  equivalent  static  load  vector  for  the  stress  response,  KA  is  the  linear 
stiffness  matrix,  aN(t)  from  Eq.  (2)  is  utilized  as  the  initial  stress  effect  in  Eq. 

(5)  for  linear  static  analysis.  z"(s)  is  the  displacement  vector  from  Eq.  (5)  and  f"(j)  is 
calculated  by  multiplying  KL  and  z %(s)  as  shown  in  Eq.  (6). 
f”  (.s')  is  used  as  follows: 

K  /, (b)z r  (.s)  =  f.”  (.s)  (7) 

The  stress  response  aL(s)  is  obtained  from  Eq.  (7)  of  linear  analysis.  Elowever,  this 
stress  response  may  not  be  exactly  the  same  as  that  from  nonlinear  analysis  because  the 
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integral  points  for  calculation  of  stresses  are  different  in  nonlinear  dynamic  analysis  and 
the  initial  stress  analysis.  The  difference  can  be  adjusted  to  aL(s)  as  follows: 

a.(s)  =  ^M  (8a) 

<?lM) 

(yjjs)  =  crjl(s)y.al(s)  (8b) 

where  a  is  the  stress  correction  factor  and  i  is  the  element  number,  a Ni(s)  is  the 
nonlinear  stress  response  from  Eq.  (2).  a Li(s)  is  the  linear  stress  response  from  Eq.  (7). 
The  stress  correction  factor  is  calculated  from  Eq.  (8a).  In  Eq.  (8b),  the  superscript  j 
means  the  iteration  number  in  linear  static  response  optimization.  The  corrected  stress 
a h(s)  is  calculated  from  Eq.  (8b).  When  j  is  equal  to  zero,  the  corrected  stress  &JLi(s) 
has  exactly  the  same  values  as  the  stress  response  from  nonlinear  dynamic  analysis  with 
the  initial  design.  The  stress  response  aJLi  (s)  is  changed  when  the  design  variables  are 
changed  in  linear  static  response  optimization.  Because  the  correction  factor  a  and  the 

equivalent  static  loads  f("  (.s)  are  constant  in  linear  static  response  optimization,  the 

corrected  stress  aJLi(s)  is  changed  as  the  design  variables  change.  Ultimately,  the 
corrected  stress  <rJLi(s)  in  linear  response  optimization  is  the  same  as  that  of  the  nonlinear 

stress  response.  Therefore,  if  the  equivalent  static  load  f "  (s)  is  used  as  an  external  load 

with  a  in  linear  static  response  optimization,  the  same  stress  as  the  one  from  nonlinear 
analysis  can  be  considered  in  the  linear  static  response  optimization  process.  Figure  3 
presents  this  process. 

If  a  problem  has  a  displacement  constraint  as  well  as  a  stress  constraint,  equivalent  loads 


should  be  calculated  with  respect  to  each  response,  and  the  sets  of  the  equivalent  static 


loads  are  utilized  in  linear  static  response  optimization  as  multiple  loading  conditions. 

2.3  The  steps  for  nonlinear  dynamic  response  structural  optimization 

using  equivalent  static  loads  (NDROESL) 

The  overall  process  of  the  NDROESL  algorithm  is  illustrated  in  Fig.  4.  The  steps  of  the 
algorithm  are  as  follows: 

Step  1.  Set  initial  design  variables  and  parameters  (design  variables:  bU)  =  b(ll),  cycle 
number:  k  =  0 ,  convergence  parameter:  a  small  number  s ). 

Step  2.  Perform  nonlinear  dynamic  analysis  with  blA  1 .  Hence  the  linear  stiffness  matrix 
and  nonlinear  responses  are  obtained. 

Step  3.  When  k  =  0,  go  to  Step  4.  When  k  >  0,  if 


to 

VI 

7 

-te 

1 

s 

-Q 

(9a) 

gj(blk+1\zN(t),aN(t))^o 

(9b) 

0  =  WO  =  !,•••,«) 


then  terminate  the  process.  Otherwise,  go  to  Step  4.  If  Eq  (9a)  is  satisfied  and 
Eq  (9b)  is  not  satisfied,  reduce  the  convergence  parameter  e  to  have  a  smaller 
value  and  go  to  Step  4. 

Step  4.  Calculate  the  equivalent  static  load  sets  as  follows: 

f;f)0)  =  Ki(b)zjV(0and  C?W0)  =  Ki(b)z"(0  (10) 

Step  5.  Solve  the  following  linear  static  response  optimization  problem: 
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Find 

b(*+1) 

(11a) 

to  minimize 

/(b<i+1)) 

(lib) 

subject  to 

Ki(b,i+1))zu)-f;’w(u=o 

(He) 

Ki(ba+1))z(U-C?W(U  =  0 

(Hd) 

s  =  1 

g/b(*+1),z(s),6(s))<0,  y=l,-,/ 

(lie) 

b|*+1)  <  b|-i+1)  <  b)y+1),  i  =  1  ,  —  ,m 

(Ilf) 

The  external  load  f  (s)  is  the  equivalent  static  load  vector  and  2 n  equivalent 
static  load  sets  are  used  as  multiple  loading  conditions  during  the  linear  static 
response  optimization  process. 

Step  6.  Update  the  design  results,  set  k  =  k  + 1  and  go  to  Step  2. 

2.4  A  small  scale  example:  nonlinear  dynamic  structural  optimization 
of  a  cantilever  plate 

A  small  scale  problem  is  solved  by  using  the  NDROESL  method  to  validate  the  method. 
The  model  is  a  cantilever  plate  with  120  shell  elements.  The  loading  and  boundary 
conditions  are  illustrated  in  Fig.  5.  Figure  6(a)  presents  the  dynamic  load  profile.  The 
duration  time  is  0.01  second  and  the  total  analysis  time  is  0.1  second.  ABAQUS  6.7  (24>  is 
used  for  the  nonlinear  dynamic  analysis.  The  implicit  method  is  used  with  a  constant 
incremental  size  of  0.0002.  The  total  number  of  time  steps  is  five  hundred  and  the  total 
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number  of  equivalent  static  loads  is  the  same.  Figure  6(b)  shows  the  strain-stress  curve  of 
the  used  material  for  this  problem.  The  material  has  bilinear  elastoplastic  strain-stress 
curve.  The  Young’s  modulus  is  68.9  GPa  and  the  tangent  modulus  is  34.5  GPa.  The 
yield  strength  is  172  MPa.  The  Poisson  ratio  is  0.35  and  the  density  is  2710  kg/m3.  Both 
geometric  and  material  nonlinearities  are  considered  in  this  problem. 

Figure  7  illustrates  the  maximum  displacement  responses  from  nonlinear  dynamic 
analysis  and  linear  static  analysis  with  ESL  for  the  initial  model.  As  shown,  the  responses 
are  exactly  the  same.  The  maximum  difference  is  4x1 0"8  m.  Therefore,  the 
transformation  is  validated.  The  optimization  formulation  is  as  follows: 


Find 

bt  0  =  1, -,29) 

(12a) 

to  minimize 

Mass 

(12b) 

subject  to 

kfp  <  20.0mm  (p  =  l,---,500) 

(12b) 

3.0  mm  <bt  <10.0  mm 

(12c) 

The  design  variables  are  the  thicknesses.  The  cantilever  plate  is  divided  into  twenty 
nine  sections  with  respect  to  the  x  direction  and  the  total  number  of  design  variables  is 
twenty  nine.  The  objective  function  is  the  mass.  The  constraint  is  that  the  magnitude  of 
maximum  displacement  should  be  less  than  the  allowable  displacement  of  20  mm  at  all  the 
time  steps. 

This  problem  is  solved  by  NDROESL  as  well  as  by  a  conventional  method.  The 
modified  method  of  feasible  directions  algorithm  in  a  commercial  optimization  code  DOT 
5.7  is  used  for  the  conventional  method.  <26)  The  finite  difference  method  (FDM)  is 
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employed  for  sensitivity  analysis.  The  results  of  both  methods  are  compared. 

Figures  8  and  9  illustrate  the  history  of  the  objective  function  and  constraint  violation  for 
NDROESL  and  FDM,  respectively.  Table  1  shows  the  optimization  results  for  the 
cantilever  plate  problem.  As  shown  in  the  table,  the  optimum  mass  is  almost  the  same. 
The  displacement  constraint  is  active  at  the  optima  of  both  methods.  Figure  10  presents 
the  maximum  displacements  from  nonlinear  dynamic  analyses  at  the  optima  of  both 
methods.  They  are  almost  the  same.  Since  the  solution  from  the  conventional  method 
can  be  considered  as  a  mathematical  optimum,  the  quality  of  the  solution  from  NDROESL 
is  excellent. 

The  efficiency  of  the  two  methods  is  quite  different.  Only  eight  nonlinear  dynamic 
analyses  are  required  in  NDROESL  while  three  hundred  and  sixty  five  analyses  are 
required  in  the  conventional  method  using  FDM.  The  same  computer,  Intel  Pentium  Dual 
CPU  3.20  GHz,  3.25  GB  RAM,  (27)  is  used  for  the  analysis  and  optimization.  In  total 
CPU  time,  NDROESL  requires  22  minutes  while  FDM  requires  486  minutes.  Figure  1 1 
illustrates  the  thickness  distribution  of  the  optimum  models  from  both  methods.  The 
thickness  of  the  root  is  thick  and  that  of  the  tip  is  thin  in  both  methods;  however,  the 
profiles  are  different.  The  difference  of  sensitivity  causes  the  difference  of  the  optimum 
profile.  The  linear  response  is  used  for  the  calculation  of  sensitivity  in  NDROESL.  On 
the  other  hand,  the  nonlinear  response  is  directly  used  for  the  calculation  of  sensitivity  in 
the  conventional  method  using  FDM.  The  difference  of  sensitivity  is  reduced  as  the  cycle 
is  repeated.  The  results  of  NDROESL  are  almost  the  same  as  those  of  the  conventional 
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method  using  FDM.  However,  NDROESL  is  more  efficient  than  the  conventional  method. 

Conceptually,  it  seems  that  the  joined-wing  structure  is  a  cantilever  type  structure.  The 
root  of  the  wing  is  fixed  at  the  fuselage.  Several  thousand  shell  elements  are  used  for  the 
finite  element  method  of  the  structure.  From  the  next  section,  the  analysis  and 
optimization  of  the  joined-wing  structure  will  be  explained.  Since  it  is  a  very  large  scale 
problem,  the  conventional  method  is  almost  impossible  to  use.  Therefore,  only  the 
NDROESL  method  is  used  for  nonlinear  dynamic  structural  optimization  of  the  joined- 
wing. 


3  Analysis  of  the  joined-wing 

3.1  Finite  element  modeling  of  the  joined-wing 

Figure  12  illustrates  a  finite  element  model  of  the  joined-wing.  The  joined-wing 
consists  of  five  parts,  which  are  the  fore-wing,  the  aft-wing,  the  mid-wing,  the  tip-wing 
and  the  edge  around  the  joined-wing.  Each  part  is  composed  of  the  top  skin,  the  bottom 
skin,  the  spar  and  the  rib.  The  length  from  the  wing-tip  to  the  wing-root  is  38  m  and  the 
length  of  the  chord  is  2.5  m.  The  model  has  3027  elements  with  2857  quadratic  elements, 
156  triangular  elements  and  14  rigid  elements.  Rigid  elements  make  connections  between 
the  nodes  of  the  aft-wing  root  and  the  center  node  of  the  aft-wing  root.  The  structure  has 
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two  kinds  of  aluminum  materials.  One  has  the  Young’s  modulus  of  72.4  GPa,  the  shear 
modulus  of  27.6  GPa  and  the  density  of  2770  kg/m3.  The  other  has  36.2  GPa,  13.8  GPa 
and  2770  kg/m3,  respectively.  The  former  material  is  used  for  the  entire  elements  except 
for  the  edge  part.  The  latter  material  is  only  used  for  the  elements  of  the  edge  part. 

3.2  Loading  conditions  of  the  joined- wing 

Eleven  static  loading  conditions  for  structural  optimization  have  been  defined  by  the 
AFRL.  (4)  These  loading  conditions  are  composed  of  seven  maneuver  loads,  two  gust 
loads,  one  take-off  load  and  one  landing  load  as  shown  in  Table  2.  Each  loading 
condition  has  a  different  loading  direction  and  magnitude.  The  gust  loading  conditions 
are  especially  important  in  these  loading  conditions.  Gust  is  the  movement  of  the  air  in 
turbulence  and  the  gust  load  has  a  large  impact  on  the  airplane.  Static  loads  for  the  gust 
can  be  generated  from  an  aeroelastic  model  which  uses  the  Panel  method.  <l2)  The  Panel 
method  is  used  to  calculate  the  velocity  distribution  along  the  surface  of  the  airfoil.  Panel 
methods  have  been  developed  to  analyze  the  flow  field  around  arbitrary  bodies  in  two  and 
three  dimensions.  The  surface  of  the  airfoil  is  divided  into  trapezoid  panels. 
Mathematically,  each  panel  generates  the  velocity  on  it.  This  velocity  can  be  expressed 
by  relatively  simple  equations  which  contain  geometric  relations,  such  as  distances  and 
angles  between  the  panels.  The  Panel  method  is  referred  to  as  the  boundary  element 
method  in  some  publications.  (12)  Detailed  explanation  of  the  Panel  method  is  out  of 
scope  of  this  work. 
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The  real  gust  load  acts  dynamically  on  the  airplane.  Also,  dynamic  loads  are  required 


for  optimization  with  equivalent  static  loads.  However,  the  generation  of  exact  dynamic 
loads  which  consider  the  nonlinear  dynamic  behavior  of  the  airplane  is  very  difficult. 
Therefore,  the  static  gust  loads  of  Reference  4  are  transformed  to  dynamic  loads. 
Generally,  there  are  several  methods  for  generating  dynamic  gust  loads.  <U)  Here,  the 
approximated  dynamic  load  is  evaluated  by  multiplying  the  static  load  by  the  1 -cosine 
function. 

The  duration  time  of  the  dynamic  gust  load  is  calculated  from  the  following  equation.  (ll) 


U  = 


U, 


de 


1  —  cos 


2  7® 

25  C 


(13) 


where  U  is  the  velocity  of  the  gust  load,  Ude  is  the  maximum  velocity  of  the  gust  load,  s 
is  the  distance  penetrated  into  the  gust  and  C  is  the  geometric  mean  chord  of  the  wing. 
The  conditions  for  the  coefficients  are  shown  in  Table  3.  From  Table  3  and  Eq.  (13),  the 
duration  time  is  0.374  seconds.  The  airplane  stays  in  the  gust  for  0.374  seconds. 

The  dynamic  gust  load  is  calculated  as  follows: 

^dynamic  =  ^static  X  f1  “  C0S  (14) 


where  F’static  is  the  static  gust  load  which  is  the  eighth  or  ninth  load  in  Table  2.  It  is  noted 
that  the  period  of  the  gust  load  is  0.374  second  and  the  duration  time  of  the  dynamic  load  is 
0.374  second.  The  dynamic  load  becomes  zero  after  0.374  second. 


3.3  Boundary  conditions  of  the  joined- wing 

The  roots  of  the  fore-wing  and  the  aft- wing  are  joined  to  the  fuselage.  That  is,  the  entire 
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part  of  the  fore-wing  root  is  attached  to  the  fuselage.  Therefore,  all  the  degrees  of 
freedom  in  six  directions  are  fixed.  On  the  other  hand,  the  aft-wing  root  can  be  rotated 
with  respect  to  the  y-axis  in  Fig.  13.  The  boundary  nodes  of  the  aft- wing  root  are  rigidly 
connected  to  the  center  node.  The  center  node  has  an  enforced  rotation  with  respect  to  the 
y-axis.  The  boundary  nodes  are  set  free  in  the  x  and  z  translational  directions.  Other 
degrees  of  freedom  are  fixed.  The  enforced  rotation  generates  torsion  on  the  aft-wing  and 
has  quite  an  important  aerodynamic  effect.  The  amounts  of  the  enforced  rotation  are  from 
-0.0897  radian  to  0.0  radian.  These  rotational  values  are  different  in  each  mission  leg. 
The  boundary  conditions  are  illustrated  in  Fig.  13. 

3.4  Nonlinear  dynamic  analysis  of  the  joined-wing 

Nonlinear  dynamic  analysis  is  performed  under  the  gust  loading  conditions.  Geometric 
nonlinearity  is  considered  in  nonlinear  dynamic  analysis.  The  dynamic  loads  are 
generated  by  Eq.  (14).  ABAQUS  6.7  (24)  is  used  for  nonlinear  dynamic  analysis.  F1P-UX 
Itanium  II  computer  is  used  for  nonlinear  dynamic  analysis <28)  As  mentioned  before,  the 
duration  time  of  the  dynamic  gust  load  is  0.374  second  and  the  total  analysis  time  is  1.8 
seconds.  The  size  of  the  time  step  is  0.1  second.  Then,  the  stress  response  is  recorded 
every  0.1  second.  Therefore,  each  loading  condition  has  18  time  steps.  Then,  the  total 
number  of  time  steps  is  thirty  six  for  the  two  gust  loading  conditions.  In  the  linear  static 
response  optimization  process  using  the  equivalent  static  loads,  thirty  six  static  loading 
conditions  are  utilized  as  multiple  loading  conditions. 
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Figure  14  illustrates  the  von  Mises  stresses  from  nonlinear  dynamic  analysis.  The  stress 
fluctuates  and  the  maximum  stress  occurs  after  0.374  second  which  is  the  duration  time  of 
the  dynamic  load.  Moreover,  the  maximum  stress  occurs  within  1.8  seconds.  Generally, 
the  maximum  stress  occurs  at  the  wing  root.  Figure  15  presents  the  stress  contour  of  the 
joined-wing  at  1.1  second  and  the  maximum  stress  occurs  under  the  gust  loading  condition 
9. 


4  Structural  optimization  of  the  joined-wing 

4.1  Definition  of  design  variables 

As  mentioned  earlier,  the  FEM  model  of  the  joined-wing  has  3027  finite  elements.  It  is 
not  reasonable  to  select  the  properties  of  all  the  elements  as  design  variables  for 
optimization.  Thus,  the  design  variable  linking  technology  is  utilized.  The  wing 
structure  is  divided  into  forty  eight  sections  and  each  section  has  the  same  thickness.  The 
finite  element  model  is  adopted  from  Reference  4.  The  joined-wing  is  divided  into  five 
parts  as  illustrated  in  Fig.  12.  Each  part  is  composed  of  the  top  skin,  the  bottom  skin,  the 
spar  and  the  rib.  Fig.  16  presents  the  division  of  the  mid-wing.  The  top  and  bottom 
skins  are  divided  into  three  sections.  The  sections  are  the  wing-skin-front,  the  wing-skin- 
middle  and  the  wing-skin-rear.  The  spars  of  the  mid-wing  are  divided  into  seven  sections. 
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The  spars  of  other  wings  are  divided  into  three  sections.  Other  parts  such  as  the  fore-wing, 
the  aft-wing,  the  wing  tip  and  the  edge  are  divided  in  the  same  manner.  The  wing  tip  and 
the  edge  parts  are  not  used  as  design  variables.  Therefore,  only  thirty  four  sections  among 
the  forty  eight  sections  are  used  as  design  variables.  Design  variables  are  defined  based 
on  Reference  9. 


4.2  Optimization  formulation 

The  optimization  problem  is  formulated  as 

Find  bt  (;  =  !,•■•, 34)  (15a) 

to  minimize  Mass  (15b) 

subject  to  \<rf  fi)|  <  <Tallowable  (15c) 

(p  =  8,9;y  =  l,---,1300;t  =  1, — ,18) 

0.001016m  <  bskinpart  <  0.05m  ( 1 5d) 

0.000254m  <Z>wingsparsandribs  <  0.08m  (15e) 


where  bt  is  the  thickness  of  the  i  th  section,  a J (t)  is  the  stress  of  the  j  th  element  at 

the  t  th  time  step  under  the  loading  condition  p .  As  mentioned  earlier,  two  gust  loading 
conditions  (the  8th  and  9th  loading  conditions  in  Table  2)  are  used  for  nonlinear  dynamic 
response  optimization. 

The  objective  function  is  the  mass.  The  mass  of  the  initial  model  is  3863  kg.  The  FE 
model  has  3027  elements  and  each  element  has  a  different  thickness.  The  design 
variables  are  linked  according  to  the  definition  of  the  design  variables.  Then  the  initial 
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mass  is  4285  kg.  The  upper  and  lower  bounds  are  defined  for  each  part.  0.001016  m 
and  0.000254  m  are  used  as  the  lower  bounds  of  the  skin  part  and  wing  spars,  and  that  of 
rib  parts,  respectively.  0.05  m  is  used  as  the  upper  bound  of  the  skin  part  and  0.08  m  is 
used  as  the  upper  bounds  of  the  spars  and  rib  parts. 

The  material  of  the  joined-wing  is  aluminum.  <4)  The  allowable  von  Mises  stress  for 
aluminum  is  set  by  269MPa.  Since  the  safety  factor  1.5  is  used,  the  allowable  stress  is 
reduced  to  179  MPa.  Stresses  of  all  the  elements  except  for  the  edge  part  and  the  wing  tip 
part  should  be  less  than  the  allowable  stress  179  MPa. 


5  Discussion 

5.1  The  results  of  nonlinear  dynamic  response  optimization 

Nonlinear  dynamic  response  structural  optimization  of  the  joined-wing  is  carried  out 
using  equivalent  static  loads.  Two  gust  loading  conditions  are  used  as  external  dynamic 
loads.  Each  loading  condition  is  divided  into  eighteen  time  steps  from  0.0  second  to  1.8 
second.  According  to  the  equivalent  static  loads  concept,  thirty  six  static  loading 
conditions  are  defined  for  the  two  gust  loads. 

Table  4  and  Fig.  17  show  the  history  of  the  optimization  process.  The  objective 
function  is  increased  by  318.4  percent  from  4285.96  kg  to  17933.55  kg.  It  is  noted  that 
the  constraints  are  satisfied  when  nonlinear  dynamic  analysis  is  performed  with  the 
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optimum  solution.  The  stress  response  of  the  optimum  is  illustrated  in  Fig.  18.  The 
critical  stresses  occur  at  0.2  second,  0.4  second,  0.9  second  and  1.4  second.  Figure  19 
illustrates  the  stress  contours  of  the  optimum  at  1.4  second.  The  maximum  stress  of 
optimum  occurs  at  element  1407  which  is  located  in  the  top  skin  of  the  aft-wing  root.  The 
magnitude  of  the  maximum  stress  is  179.9  MPa  at  the  time  of  1.4  second.  Generally,  the 
effect  of  loading  condition  9  (cruise  speed  gust  load)  is  more  severe  than  that  of  loading 
condition  8  (maneuver  speed  gust). 

The  violation  of  the  stress  constraint  of  cycle  7  is  smaller  than  that  of  cycle  8.  Flowever, 
the  mass  of  cycle  8  is  smaller  than  that  of  cycle  7.  The  process  is  considered  as 
converged  when  the  difference  between  the  design  variables  of  the  current  cycle  and  those 
of  the  previous  cycle  is  smaller  than  a  given  small  number.  The  convergence  criteria  are 
satisfied  in  cycle  8.  Both  results  of  cycle  7  and  8  may  be  selected  as  the  optimum  design. 

5.2  Discussion  about  the  optimum  design 

As  mentioned  earlier,  the  mass  is  increased  by  318.4  percent  from  4285.96  kg  to 
17933.55  kg.  Overall,  the  optimum  thickness  from  nonlinear  dynamic  response  structural 
optimization  is  larger  than  that  of  the  initial  model.  In  Reference  9,  where  linear  dynamic 
response  optimization  of  a  joined-wing  is  performed,  the  optimum  mass  is  12725.52  kg. 
The  optimum  mass  of  nonlinear  dynamic  response  optimization  is  larger  than  that  of  linear 
dynamic  response  optimization.  It  is  reasonable  because  the  geometric  nonlinearity  is 
added  in  this  research.  The  stress  constraint  violation  of  the  initial  model  is  344.21%  in 
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Reference  9.  However,  the  stress  constraint  violation  of  the  initial  model  is  736.9%  in 


this  research.  The  mass  of  the  initial  design  and  the  definition  of  the  design  variables  are 
not  exactly  the  same  between  Reference  9  and  this  research.  However,  the  FEM  model, 
the  boundary  conditions  and  the  critical  loading  conditions  are  the  same.  Therefore,  it 
seems  that  this  comparison  is  useful  for  the  design  of  a  joined-wing. 

Table  5  and  Fig.  20  show  the  optimum  thickness.  The  thicknesses  of  the  parts  in  the  aft- 
wing  are  quite  large.  The  leading  edge  and  middle  part  of  the  top  skin,  the  trailing  edge 
part  of  the  bottom  skin  and  the  leading  and  trailing  spars  of  the  aft-wing  are  as  thick  as  the 
upper  bounds.  In  Fig.  20,  dv/  means  the  design  variable  number.  In  the  top  skin  of  the 
aft-wing,  the  thicknesses  of  the  leading  edge  and  middle  section  are  larger  than  that  of  the 
trailing  edge.  On  the  other  hand,  in  the  bottom  skin  of  the  aft-wing,  the  thickness  of  the 
trailing  edge  is  larger  than  that  of  the  leading  edge.  This  means  that  a  large  torsion  force 
occurs  at  the  aft-wing  root  position.  This  torsion  effect  is  observed  in  Fig.  19.  In  the  top 
skin  of  the  aft-wing  root,  a  large  stress  occurs  at  the  leading  edge.  At  the  bottom  skin  of 
the  aft-wing  root,  a  large  stress  occurs  at  the  trailing  edge.  It  seems  that  the  optimum 
thicknesses  of  the  aft-wing  are  influenced  by  the  torsion  effect. 


6  Conclusions 
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The  joined-wing  is  defined  as  an  airplane  that  incorporates  tandem  wings  arranged  to 
form  diamond  shapes  in  both  top  and  front  views.  The  joined-wing  configuration  has 
many  advantages  from  the  viewpoint  of  aerodynamic  performance  and  weight  reduction. 
However,  due  to  the  specific  shape  of  the  joined-wing,  it  has  large  geometric  nonlinearity 
under  the  gust  loading  conditions.  The  real  gust  acts  dynamically.  Therefore,  the 
nonlinear  dynamic  behavior  should  be  considered  in  structural  optimization  of  a  joined- 
wing.  The  dynamic  gust  load  profile  is  calculated  by  multiplying  the  static  gust  loads  by 
the  1 -cosine  function. 

The  equivalent  static  loads  are  used  for  nonlinear  dynamic  response  structural 
optimization  of  a  joined-wing.  The  existing  concepts  of  the  equivalent  static  loads  are 
expanded  for  nonlinear  dynamic  response  optimization.  This  is  called  nonlinear  dynamic 
response  optimization  using  equivalent  static  loads.  (NDROESL)  The  equivalent  static 
loads  are  defined  as  the  linear  static  load  sets  which  generate  the  same  response  field  in 
linear  static  analysis  as  that  from  nonlinear  dynamic  analysis.  Therefore,  if  equivalent 
static  loads  are  used  as  applied  loads,  the  same  responses  from  nonlinear  dynamic  analysis 
can  be  considered  in  linear  static  response  optimization.  Equivalent  static  loads  are  made 
to  reduce  the  number  of  nonlinear  dynamic  analyses.  Also,  because  the  method  is 
gradient  based  optimization,  the  solution  is  exact.  An  example  of  the  cantilever  plate  is 
solved  by  NDROESL  as  well  as  the  conventional  method  using  the  finite  difference 
method.  By  comparing  results,  NDROESL  is  more  efficient  than  the  conventional 
method.  And  the  objective  function  values  of  the  two  methods  are  almost  the  same. 
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The  optimum  design  considering  the  nonlinear  dynamic  effect  of  the  joined-wing 
satisfies  all  the  stress  constraints.  The  joined-wing  is  divided  into  forty  eight  sections  and 
the  thicknesses  of  thirty  four  sections  are  used  as  design  variables  for  nonlinear  dynamic 
response  optimization.  The  mass  is  increased  by  318.4  percent.  It  is  because  the 
constraint  violation  of  the  initial  model  is  quite  large  and  the  thicknesses  of  almost  all  the 
sections  are  increased.  It  is  noted  that  the  optimum  thicknesses  of  the  leading  edge  of  the 
top  skin  and  those  of  the  trailing  edge  of  the  bottom  skin  in  the  aft- wing  are  5  cm,  which  is 
the  upper  bound.  Only  eight  nonlinear  dynamic  analyses  are  required  for  nonlinear 
dynamic  response  structural  optimization  of  the  joined-wing.  Nonlinear  dynamic 
response  optimization  of  a  joined-wing  using  the  proposed  method  is  very  successful  and 
efficient  although  the  problem  is  fairly  large  scale. 
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Fig.  1  Configuration  of  the  joined- wing 
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Fig.  2  Generation  of  equivalent  static  loads  for  displacement  constraints 
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Fig.  3  Generation  of  equivalent  static  loads  for  stress  constraints 
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Fig.  4  Optimization  process  using  the  equivalent  static  loads 
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Fig.  5  Cantilever  plate  structure 
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(a)  Dynamic  load  profile  of  the  cantilever  plate  (b)  Bilinear  elastoplastic  strain-stress  curve 
structure 


Fig.  6  Dynamic  load  profile  and  strain-stress  curve  of  the  cantilever  plate  structure. 
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Fig.  7  Displacement  response  from  nonlinear  dynamic  analysis  and  linear  static  analysis 

with  ESL  for  the  initial  model 
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Fig.  8  Objective  function  and  constraint  violation  history  of  the  cantilever  plate  structure 

using  NDROESL 
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Fig.  9  Objective  function  and  constraint  violation  history  of  the  cantilever  plate  structure 

using  FDM 
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Objective  value  —♦—Constraint  violation 


Fig.  1 0  Displacement  response  from  nonlinear  dynamic  analysis  with  optimum  design  of 

NDROESL  and  FDM 
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Fig.  11  Thickness  distribution  of  the  optimum  models  from  NDROESL  and  FDM 
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Tip- wing 


Fig.  12  Finite  element  modeling  of  the  joined-wing 
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Fig.  13  Boundary  conditions  of  the  joined-wing 
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Fig.  14  Stress  response  of  the  joined-wing  under  the  cruise  speed  gust  loading  condition 


(a)  Top  skin 

Fig.  15  Stress  contours  at 


(b)  Rib  and  spar  (c)  Bottom  skin 

second  from  nonlinear  dynamic  analysis  of  a  joined-wing 
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Fig.  16  Sections  for  definition  of  design  variables 
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Fig.  17  Flistory  of  nonlinear  dynamic  response  structural  optimization 
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Fig.  18  Stress  response  of  the  optimum  design  under  the  cruise  speed  gust  loading 

condition 


Fig.  19  Stress  contours  of  the  optimum  design  at  1.4  second 
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Fig.  20  Thickness  contour  of  the  nonlinear  dynamic 


(c)  Bottom  skin 
response  optimization  result 
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Table  1  Optimum  results  for  the  cantilever  plate  problem 


Initial 

Optimum  of 

NDROESL 

Optimum  of 

FDM 

Mass 

0.09756  kg 

0.13353  kg 

0.13314kg 

Maximum  displacement 

59.6  mm 

20.00  mm 

20.06  mm 

Number  of  iterations  (cycles) 

8 

9 

Number  of  nonlinear  transient 

analyses 

8 

365 

Number  of  nonlinear  transient 

analyses  except  for  gradient  call 

104 

Total  number  of  iterations  for  linear 

response  optimization 

24 

Total  CPU  time 

22  minutes 

486  minutes 
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Table  2  Load  data  of  the  joined-wing 


Number  of  loading  condition 

Load  type 

Mission  leg 

1 

2.5  g  PullUp 

Ingress 

2 

2.5  g  PullUp 

Ingress 

3 

2.5  g  PullUp 

Loiter 

4 

2.5  g  PullUp 

Loiter 

5 

2.5  g  PullUp 

Egress 

6 

2.5  g  PullUp 

Egress 

7 

2.5  g  PullUp 

Egress 

8 

Gust  (Maneuver) 

Descent 

9 

Gust  (Cruise) 

Descent 

10 

Taxi  (1.75  g  impact) 

Take-off 

11 

Impact  (3.0  g  landing) 

Landing 

Table  3  Aerodynamic  data  for  the  joined-wing 


Gust  maximum  velocity 

18.2  m/s 

Flight  velocity 

167  m/s 

Geometric  mean  chord  of  wing 

2.5  m 

Distance  penetrated  into  gust 

62.5  m 
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Table  4  Results  of  nonlinear  dynamic  response  optimization  of  the  joined-wing 


Iteration  no. 

Optimum  value  (kg) 

Constraint  violation  (%) 

0 

4285.96 

736.9 

1 

10078.70 

54.2 

2 

13608.93 

19.8 

3 

15857.66 

15.0 

4 

17129.80 

8.9 

5 

17326.98 

6.7 

6 

17797.18 

2.3 

7 

17944.14 

0.0 

8 

17933.55 

0.5 
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Table  5  Optimum  thicknesses  from  nonlinear  response  optimization  using  ESL 


Number 

Name 

Initial  thickness 

(meter) 

Optimum  thickness 

(meter) 

1 

AFTBOTTOMSKINFRONT 

0.001020 

0.002626 

2 

AFTBOTTOMSKINMIDDLE 

0.001554 

0.040679 

3 

AFTBOTTOMSKINREAR 

0.023830 

0.049995 

4 

AFTRIB 

0.000298 

0.010910 

5 

AFTSPARFRONT 

0.000288 

0.079995 

6 

AFTSPARMIDDLE 

0.000798 

0.048723 

7 

AFTSPARREAR 

0.080000 

0.079995 

8 

AFTTOPSKINFRONT 

0.041890 

0.049995 

9 

AFTTOPSKINMIDDLE 

0.002718 

0.049996 

10 

AFTTOPSKINREAR 

0.001493 

0.017102 

11 

FOREBOTTOMSKINFRONT 

0.001442 

0.005083 

12 

FOREBOTTOMSKINMIDDLE 

0.003319 

0.006304 

13 

FOREBOTTOMSKINREAR 

0.009598 

0.013159 

14 

FORERIB 

0.000274 

0.000356 

15 

FORESPARFRONT 

0.000615 

0.000479 

16 

FORESPARMIDDLE 

0.000732 

0.001265 

17 

FORESPARREAR 

0.012790 

0.009565 

18 

FORE  TOP  SKIN  FRONT 

0.002578 

0.003251 

19 

FORE  TOP  SKIN  MIDDLE 

0.004656 

0.009070 

20 

FORETOPSKINREAR 

0.001539 

0.003532 

21 

MIDBOTTOMSKINFRONT 

0.001020 

0.004202 

22 

MIDBOTTOMSKINMIDDLE 

0.002000 

0.003295 

23 

MIDBOTTOMSKINREAR 

0.001261 

0.002538 

24 

MIDRIB 

0.000384 

0.001997 

25 

MIDSPAR 

0.000838 

0.000398 

26 

MIDSPARl 

0.000254 

0.001663 

27 

MIDSPAR2 

0.000254 

0.000254 

28 

MIDSPAR3 

0.000254 

0.000254 

29 

MIDSPAR4 

0.000542 

0.003763 

30 

MIDSPAR5 

0.000254 

0.000339 

31 

MIDSPAR6 

0.000458 

0.016899 

32 

MIDTOPSKINFRONT 

0.001020 

0.005426 

33 

MIDTOPSKINMIDDLE 

0.001793 

0.004075 

34 

MIDTOPSKINREAR 

0.001020 

0.001090 
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